Get the top 10 brand names

top_10_brands <- stations_2019 %>% 
  count(brand) %>%
  arrange(desc(n)) %>% 
  filter(!is.na(brand)) %>% 
  head(10) %>% 
  pull(brand)

top_10_brands

##  [1] "ARAL"       "Shell"      "ESSO"       "TOTAL"      "AVIA"      
##  [6] "JET"        "STAR"       "Agip"       "Raiffeisen" "HEM"

# Filter stations and transform longitude and latitude
geo <- stations_2019 %>% 
  filter(brand %in% top_10_brands[1:5]) %>%
  filter(latitude > 0 & longitude > 0) %>% 
  sample_frac(0.1) %>% 
  sf::st_as_sf(coords = c("longitude", "latitude"), crs=4326) %>% 
  sf::st_transform(crs=4326)

# Combine data with coordinates
station_coords <- cbind(geo, sf::st_coordinates(geo))

# Plot german map and points from station_cords on top
plt <- ggplot() + 
  geom_sf(data=landkreise_in_germany, alpha=0.2, colour="grey") +
  geom_point(data=station_coords, aes(X, Y, fill=brand), alpha=0.8) +
  theme_bw()

plotly::ggplotly(plt)